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The BABAYAGA event generator 

CM. Carloni Calame a b , G. Montagna b a , 0. Nicrosini a b and F. Piccinini a b . 

a Istituto Nazionale di Fisica Nucleare, Sezione di Pavia, Via A. Bassi, 6 - 27100 Pavia, Italy 
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The program BABAYAGA is an event generator for QED processes at flavour factories, mainly intended for luminosity 
measurement of e + e~ colliders in the center of mass range 1-10 GeV. Recently, the tt + -k~ channel has been added 
as well. The relevant (photonic) radiative corrections are simulated by means of a Parton Shower in QED. The 
theoretical precision of the approach is estimated and some phenomenological results are discussed. 



1. Introduction 

The precise determination of the machine lumi- 
nosity is an important ingredient for the success- 
ful achievement of the physics programme at the 
e + e~ colliders running with center of mass energy 
in the range of the low lying hadronic resonances 
(1 - 10 GeV). 

One of the most important challenges is the 
precise measure of the R ratio, by means of the 
energy scan or the radiative return method. The 
aim is to reduce the theoretical error on the 
hadronic contribution to the vacuum polariza- 
tion, which will reflect on the error of the anoma- 
lous magnetic moment of the muon and the 
QED coupling constant at the Z peak o.qed(Mz) 
PQ. Actually, the R measurement is going to be- 
come a precision measurement |2IH] and it will 
give a stringent test of the Standard Model pre- 
dictions. 

In this perspective, a precise knowledge of the 
cross section for QED processes for luminometry, 
expecially for Bhabha scattering, and their Monte 
Carlo (MC) simulation are crucial. The theoreti- 
cal accuracy should be better than 0.5%. 

2. The theoretical approach in BABAYAGA 

BABAYAGA 4 5 is a Monte Carlo event gener- 
ator for the simulation of QED and Tr + n~ pro- 
cesses at flavour factories. 

In order to achieve an accurate prediction of the 
cross sections of the QED processes, the relevant 



radiative correction (RC) have to be accounted 
for, and, in particular, photon radiation effects 
must be included in the calculation. 

2.1. The QED Structure Functions 

The corrected cross section and the event gen- 
eration can be obtained according to the following 
master formula (see 0] and references therein): 

<r(s) = J dx\dx 2 dy\dy 2 J dfl x 

D(x 1 ,Q 2 )D(x 2 ,Q 2 )D(y 1 ,Q 2 )D(y 2 , Q 2 ) x 
da (xix 2 s) 



dQ 



-Q(cuts) 



(1) 



The previous equation is used in BABAYAGA for the 



generation of the processes e H 



("7), 



— > n + n~{n^), — > 77(717) and — * 7r + ir~ (wy) . The 
last channel has been recently implemented and, 
in this case, only photons coming from initial 
state are accounted for. 

In eq. n the Born-like cross section for the 
process dao/dft is convoluted with the QED 
non-sing let Structure Functions (SF's) D(x,Q 2 ), 
which are the solution of the Dokshitzer-Gribov- 
Lipatov-Altarelli-Parisi (DGLAP) [H] equation in 
QED 

Q 2 lL D M 2 ) = ^ f 1 -P + (y)D(-,Q 

0Q 1 2ir J x y y 



P + (x) = 
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-«y(l-x) / dtP(t) (2) 



1 - x jo 
where P+(x) is the regularized Altarelli-Parisi 
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vertex. The SF's in eq. 2] account for photon 
radiation emitted by both initial-state (ISR) and 
final-state (FSR) fermions. The SF's allow to in- 
clude the universal virtual and real photonic cor- 
rections, up to all orders of perturbation theory. 
Strictly speaking, radiation in eq. ^ is collinear 
to the emitting fermions. 

It is possible to go beyond the collinear approxi- 
mation by solving the DGLAP equation by means 
of a MC algorithm, the so-called Parton Shower 
algorithm (PS), which is discussed in the next 
section. 

2.2. The Parton Shower algorithm 

The Parton Shower is a MC algorithm which 
solves exactly the DGLAP equation (see refer- 
ences in 0]) and it originates from an iterative 
solution of eq. 



D(x, Q 2 ) = II(Q 2 , m 2 )5(l - x) 
ct [ s „, , s ds' 

27T 7^2 



+— / n(QV)^ri(s',TO 2 ) x 



x / dyP(y)5(x-y) + 
Jo 

+ 2 photons + • • • 



(3) 



In eq. □ II(Q 2 ,m 2 ) = log ^ £ p ^ is the 
Sudakov Form Factor, where the soft and virtual 
RC are exponentiated. The Q 2 scale is not a pri- 
ori fixed, but it has to be chosen in order to repro- 
duce as close as possible the leading RC contribu- 
tions coming from the exact order a calculation. 
For example, for the Bhabha scattering, a reason- 
able choice is Q 2 = st/ue, where s,t and u are the 
usual Mandelstam variables and e is the Euler's 
number. 

The main advantage of this solution is that it 
keeps track of the number of the emitted photons 
and that their kinematics can be approximately 
recovered. Therefore, the whole branching event 
can be generated and the multi-photon emission 
can be exclusively simulated. 

The generation of the photon's momenta has 
to be carefully considered. As suggested by eq. 
13 when the PS generates photons, their energy 
is extracted from the Altarelli-Parisi vertex P(x), 
where x is the energy fraction remaining to the 
fcrmion after the photon emission. 



In order to generate also the angular variables, 
different solutions can be exploited [S], being 
those degrees of freedom integrated out in eqs. 
121 and |3J The first possible choice is to weight 
the photon angles according to the leading poles 
which are present in the radiation matrix element. 
Namely, the photons' angles can be distributed 
according to the following formula 



cos 



E 



i 



pi 



k 



(4) 



where pi are the momenta of the charged particles 
and k is the momentum of the generic photon. If 
eq. 0] can be a first approximation, it is not so 
accurate when photonic variables are looked at 
closer. In order to improve the exclusive photon 
generation, we can observe that in eq. 01 the in- 
terference effects between ISR and FSR are com- 
pletely neglected, thus the fermions radiate inde- 
pendently. 

The radiation coherence effects can be included 
if we consider the cross section for the produc- 
tion of n photons from a generic process involv- 
ing fermions as external particles. Taking the soft 
limit and keeping all the interference terms, the 
resulting formula is 



mvjPi -Pj 

i —7 (Pi ■ k)(pj ■ h) 



(5) 



where fc; and pj are photons and fermions mo- 
menta and rji = ±1 is a charge factor. The an- 
gular spectrum of the photons can be directly ex- 
tracted from the eq. El This highly enhances 
the reliability of the PS also in the simulation of 
radiative events, as will be shown in section 0] 

An important point is the estimate of the phys- 
ical precision of the PS approach. In order to give 
a consistent estimate, an 0(a) PS has been devel- 
oped as well, which calculates the corrected cross 
section of eq. ^ at first order in a in the PS ap- 
proximation. By a systematic comparison of the 
the 0(a) PS results to the exact 0(a) matrix ele- 
ment ones, the PS precision can be estimated and, 
in particular, the size of the O(a) non- leading-log 
contributions missing in a pure PS approach can 
be quantified. 
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2.3. The hard scattering cross sections 

The other ingredient entering eq. ^is the "ker- 
nel" cross section, which is the Born-like cross 
section for the process evaluated at the reduced 
center of mass energy x\x?,s. In BABAYAGA, the 
standard formulae for Bhabha, 77 and pro- 
cesses are implemented. For the 7r + 7r _ channel, 
the pion Form Factor as parameterized in |Jj is 
present. The hard scattering cross section in- 
cludes also vacumm polarization contributions. 
It is worth discussing in some details how they 
are treated in the last official release of BABAYAGA 
(3.5). 

The running of a(q ) is given by 



a(q 2 ) = 



a(0) 



1 - Aa(q 2 



Aatf) = Aa^ d tf) + Aa t tf 

+ YA ai (q 2 ) 



(6) 



The lepton and top quark contributions (Aai(q 2 ) 
and Aat(q 2 )) do not show any problem, being 
analytically calculable for any q 2 value. 

The hadronic contribution Aa^ d (q 2 ) is calcu- 
lated by means of the subroutine HADR5 8 (which 
is now superseded by a new version). This sub- 
routine is not suited to work in the q 2 time-like 
region below (40 GeV) 2 . As suggested by one of 
the authors the following solution was there- 
fore adopted: in this q 2 region, the q 2 is flipped 
to — q 2 , in order to include at least the leading 
terms of the hadronic contribution to Aa{q 2 ). 

Recently, a new version of the routine has been 
written with an improved treatment of the 
resonances region (including new experimental 
data) and which is suited for any q 2 value. It has 
been preliminarly included in the program, show- 
ing a difference on the integrated Bhabha cross 
section at DA$NE energies at a few per mille 
level. However, a deeper look is needed and will 
be carried out in the next future. 

3. The benchmark calculation LABSPV 

One of the crucial points to be considered is 
the estimate of the theoretical accuracy of the PS 
approach. In order to quantify the theoretical er- 



ror, the independent program LABSPV has been 
derived from the code SABSPV JT] , which was ex- 
tensively used for LEP luminometry. It is not 
a true event generator, but it is a cross section 
integrator. 

In LABSPV, the exact order a matrix element 
is used, and the C(a)-corrected cross section is 
calculated according to the standard formula 



(J^ +V {E 1 <k Q ) + a%{E 7 >k ) 



(7) 



where (7g +v is the soft-photon plus virtual cross 
section, a^j is the hard-photon emission cross sec- 
tion and fco is a soft-hard radiation separator, 
which the total cross section does not depend on. 
The code LABSPV calculates the cross sections for 
the Bhabha, 77 and (only for high energies, 

E c .m. ^> m M ) processes, using matrix elements 
and formulae which can be found in [12113114] . 
In the Bhabha and /i + pT case, vacuum polariza- 
tion is treated as done in BABAYAGA. 

On top of the exact 0(a) cross section of eq. 
in LABSPV, higher order leading-log QED cor- 
rections in a SF approach can be also included, 
according to the relation 



h.o. 



a SF + a SF 



(8) 



where <jg F is the order a expansion of eq. ^and 
a'gp is eq. ^ In LABSPV, the SF's are an analytical 
solution of the QED DGLAP equation ^Hj and 
are strictly collinear. Equation |H1 can be further 
improved if cast in the factorized form 



h.o. 



«1 



'SF - 



O"0 



(9) 



It can be shown JSj that in expression[5]the bulk 
of the a 2 L corrections (L is the collinear loga- 
rithm) are also taken into account, in addition 
to the exact order a and leading-log higher order 
contributions. As a consequence, the expected 
theoretical accuracy of this approach for cross sec- 
tion calculation can be estimated at the level of 
0.1%. 

Being used as a "benchmark" , LABSPV has been 
compared to different independent calculations 
finding very good agreement. For example, in 
table m the total cross section for 77 is compared 
to the analytical results of |13| . 
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Table 1 

Comparison at different center of mass energies 
of the total 77 cross section between LABSPV and 
the results of ref. [T3] (B. & K.) 



e + e — ► 77(7) total cross sections 



energy (GeV) 


B. & K. 


LABSPV 


1 


2013.7 


2013.6(3) 


3 


264.90 


264.92(5) 


5 


102.58 


102.59(2) 


7 


54.82 


54.82(1) 


10 


28.184 


28.182(6) 



4. Phenomenological results 

4.1. Coherent Parton Shower 

As discussed in paragraph 12. 21 the photon an- 
gular generation can be performed in order to 
simulate also the ISR-FSR interference, which 
is neglected in a pure "leading-pole" generation. 
The coherence effect can be important when ex- 
clusive photon distributions are considered. 

In fig. ^ the photon angular distribution is 
shown for e + e~ — ► e + e~7 events at y/s — 1.02 
GeV, with 20° < 0± < 160°, E± > 0.4 GeV 
and E 1 > 10 MeV. The leading-pole and coherent 
0{a) PS results obtained with BABAYAGA are com- 
pared to the exact 0(a) distribution of LABSPV. 
The improvement of the coherence PS is evident. 
In fig. [3 the distribution of the photon energy 
is plotted for the same sample of events. Even if 
the leading-pole PS already produces a good dis- 
tribution, being the photon energy driven by the 
Altarclli-Parisi vertex, the coherent PS is better 
also in the tail of the energy spectrum. 

A similar behaviour is found for the PS at all 
orders: in fig. the angular distribution of the 
most energetic photon as obtained by the coher- 
ent and the leading-pole PS at all orders are com- 
pared to the exact 0(a) distribution (even if this 
is not fully consistent). 

4.2. Theoretical error 

The theoretical error of BABAYAGA for cross sec- 
tion calculation can be evaluated by comparing 




20 40 60 80 "00 120 140 160 180 

photon angle (°) 



Figure 1. The photon angular distribution as ob- 
tained by the leading-pole 0(a) PS (red dots), 
the coherent 0(a) PS (blue dots) and the exact 
O(a) matrix element (histogram). Cuts are given 
in the text. 



the 0(a) PS to the exact matrix element of eq. 
{7\ Actually, the main contribution to the error 
comes from the non-leading-log 0{a) terms, miss- 
ing in the PS approach. 

In fig. the relative deviation between the 
inclusive Bhabha cross section obtained with 
BABAYAGA at 0(a) and the exact matrix element 
is plotted as a function of the acollinearity cut, 
both for the leading-pole PS and the coherent PS. 
The imposed experimental cuts are s/s — 1.019 
GeV, 20° < i?± < 160° (top in the figure) or 
50° < i?± < 130° (bottom in the figure), E^ ln = 
0.4 GeV. 

From fig. 0] it emerges that the error in 
BABAYAGA is at level of 0.5% or lower. It is worth 
noticing that for the coherent PS the size of the 
missing O(a) is almost constant as the cuts vary, 
while for the leading-pole PS it is not constant 
and it can become larger than 0.5%. 
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Figure 2. The photon energy distribution as ob- 
tained by the leading-pole 0{a) PS (red dots), 
the coherent 0(a) PS (blue dots) and the exact 
O(a) matrix element (histogram). Cuts are given 
in the text. 



Figure 3. The most energetic photon distribution 
as obtained by the leading-pole PS (open circles) 
and the coherent PS at all orders (stars) com- 
pared exact to the O(a) matrix element photon 
distribution (histogram) . 



For the Bhabha case, we can fix the theoretical 
error of BABAYAGA in the calculation of the QED 
corrected cross section at the level of 0.5%. 

4.3. Higher order corrections 

Aiming at a high precision calculation of the 
QED processes cross section, higher order con- 
tributions can not be neglected, expecially if a 
"bare" event selection criterium is adopted and 
almost-elastic events are selected. 

The PS allows to include higher-order correc- 
tions in a natural way. In fig. the relative dif- 
ference between the cross sections obtained with 
BABAYAGA at 0{a) and at all orders is shown. 
The simulation has been performed applying the 
same cuts as in fig. 0J As it can be seen, the 
higher-order contributions are quite large: they 
can change the cross section for an amount of 
1 — 2%, depending on the cuts, and are there- 



fore important to achieve the desired theoretical 
accuracy. The higher-order impact can be also es- 
timated with LABSPV, finding a good agreement. 
In the plots of figure the effect of the a 2 L con- 
tributions is also shown, estimated with LABSPV. 
They are negligible, being below the 5 • 10 -4 . 

It interesting to notice that multi-photon emis- 
sion has a sizeable effect also on distributions. For 
example, in figure H3 the distribution of radiation 
energy, defined as E c m — E + — E_ , is plotted and 
the O(a) PS is compared to the all order PS. 

5. Other channels 

BABAYAGA can generate also /i + fi~ (717) , 77(717) 
and 7r + 7r~ (717) final states, generating events ac- 
cording to eq. and the PS algorithm. 

For /tt + /i _ and 77 processes, a study as detailed 
as for Bhabha process to establish the theoretical 



Figure 4. Relative effects of missing 0(a) con- 
tributions in BABAYAGA. Blue triangles represent 
the coherent PS and the red ones the leading-pole 
PS. See text for details. 



Figure 5. Relative effects of higher-order contri- 
butions obtained with BABAYAGA and a 2 L contri- 
bution obtained with LABSPV. See text for details. 



accuracy has not been performed yet, but it is 
planned. 

As a general consideration, it should be noticed 
that in the benchmark program LABSPV the fully 
massive /^ + /i~(7) exact 0{a) matrix element is 
not present and therefore it is not reliable at low 
energies, where the muon mass effects become im- 
portant. 

Concerning the 77 final state, it has to be 
pointed out that the PS approach suffers for dou- 
ble counting. In fact, the photons generated by 
the PS are independent from the photons of the 
hard scattering process: if the PS photons are al- 
lowed to go inside the phase space region of the 
hard scattering ones, this leads to a double count- 
ing error. However, the error can be kept small 
if large angle photons are selected, because PS 
emits photons preferably collinear to initial-state 
e + and e~ . 

In the version 3.5 of the program, also the 



7T + 7r (717) channel has been added. In this case, 
only ISR is simulated with the PS. With the PS 
approach, also the radiative 7r + 7r~7 event can be 
studied, which is extremely important for the R 
measurement via the radiative return method. 
The process was added only for cross-checking 
purposes, being the PHQKHARA MC event gener- 
ator |17| much more accurate for 7r + n~"f. 

Nevertheless, the results of BABAYAGA are quite 
good. In figureEl the tt + tt~ invariant mass distri- 
bution obtained by BABAYAGA is compared to the 
PHQKARA 1.0 distribution. The imposed cuts are 
y/s = 1.02 GeV, 55° < tf± < 125°, 5° < i? 7 < 21° 
and E 7 > 10 MeV. 

6. Towards a new release of BABAYAGA 

As discussed in the section the main source 
of error in a PS approach comes from the miss- 
ing O(ot) non-log terms. It would be extremely 
useful to include them in the generator, in or- 
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Table 2 

Preliminary comparisons for Bhabha cross sections (in nb) obtained with LABSPV, BABAYAGA 3.5 and an 
improved version of BABAYAGA. Cuts are given in the text. 



acceptance 


RC 


LABSPV (eq. |HJ| 


LABSPV (eq. EJ) 


BABAYAGA 3.5 


improved BABAYAGA 


55° - 125° 


0(a) 


451.70(1) 




454.76(1) 


451.75(1) 




C(a)+h.o. 


456.20(2) 


456.36(4) 


458.86(1) 


455.96(2) 


20° - 160° 


0(a) 


6061.9(2) 




6088.8(1) 


6062.5(2) 




O(a)+h.o. 


6087.1(3) 


6090.0(7) 


6114.4(1) 


6088.2(3) 




pion system Q 2 (GeV 2 ) 



Figure 6. Effect of multi-photon emission on the 
radiation energy distribution 



Figure 7. ir + ir invariant mass distribution by 
BABAYAGA and PH0KHARA 1.0. Cuts are given in 
the text. 



der to further reduce its theoretical error. In the 
measurements of the KLOE collaboration for the 
R ratio, the theoretical error in the luminosity is 
going to be the dominant one [3] . 

When trying to include the exact O(a) matrix 
element in the PS for exclusive event generation 1 , 
different merely technical problems arise (nega- 

1 This is an extensively studied topic also in QCD. See for 
instance 1181 . 



tive weigths, higher-order not well under control, 
etc.). Work in this direction is in progress. 

Here, we would like to show some preliminary 
results obtained with a new version of BABAYAGA, 
where part of the exact 0{a) matrix element for 
Bhabha scattering is used. 

Roughly, the idea is the following: the events 
are reweigthed for a (infrared-independent) factor 
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which accounts for the difference between the ex- 
act and the PS soft-plus- virtual cross section. Ra- 
diative events should then be reweigthed with the 
exact 0(a) radiative matrix element. At present, 
we did not succeed to perform exactly the latter 
step. Nevertheless, it is possible to find out from 
eq. a weight factor which is much closer to the 
exact matrix element than what is implemented 
in BABAYAGA 3.5. The technical details will be 
discussed elsewhere. 

It is worth stressing that this procedure is done 
at differential level and the order a expansion of 
the resulting cross section reproduces exactly the 
soft-plus- virtual cross section (<jg +v ) of eq. Hand 
"almost" exactly the hard photon term (o~ r £ I ) of 
eq. □ 

In table El the cross section of LABSPV, 
BABAYAGA 3.5 and new version are compared. The 
applied cuts for the cross sections shown in the 
table are y/s = 1.02 GeV, 20° < ?9± < 160° or 
55° < 0± < 125°, E£ in > 0.408 GeV, maximum 
acollinearity £ max = 10°. In tabled the vacuum 
polarization contributions have been switched off. 
The differences between the improved version of 
BABAYAGA and LABSPV are below the 0.1%. 



7. Conclusions 

The program BABAYAGA is a MC event gener- 
ator for QED processes and 7r + 7r~ final state at 
flavour factories. The QED radiative corrections 
are included by means of a Parton Shower. The 
theoretical accuracy of the approach is estimated 
to be at 0.5% for cross section calculation. 

Some preliminary results have been presented 
of a new version of the generator, which is under 
development. In the next future, it will include 
also the exact 0(a) matrix element and will im- 
prove the theoretical accuracy of the approach. 
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